home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
Cream of the Crop 1
/
Cream of the Crop 1.iso
/
PROGRAM
/
SNIP1292.ARJ
/
RAYS.C
< prev
next >
Wrap
C/C++ Source or Header
|
1991-09-23
|
5KB
|
180 lines
/*
** public domain ray tracing functions by Daniel Lyke
*/
#include <math.h>
#include "plot.h"
struct RAY
{
double dx, dy, dz; /* Direction vector */
double ox, oy, oz; /* Origin */
};
struct PLANE
{
double nx, ny, nz; /* Vector normal (perpendicular) to plane */
double px, py, pz; /* Point on plane */
};
struct SPHERE
{
double cx, cy, cz; /* Center of sphere */
double r2; /* Radius squared */
};
struct VECTOR
{
double dx, dy, dz; /* Three dimensional vector */
};
/*
** Return the closest point of intersection of a ray with a sphere.
** Rays starting within the sphere are returned as not intersecting.
*/
double sphere_intersect(struct RAY ray, struct SPHERE sphere)
{
double a, b, c, t1, t2, t3, close, farther;
a = ray.dx * ray.dx + ray.dy * ray.dy + ray.dz * ray.dz;
close = farther = -1.0;
if(a)
{
b = 2.0 * ((ray.ox - sphere.cx) * ray.dx
+ (ray.oy - sphere.cy) * ray.dy
+ (ray.oz - sphere.cz) * ray.dz);
c = (ray.ox - sphere.cx) * (ray.ox - sphere.cx)
+ (ray.oy - sphere.cy) * (ray.oy - sphere.cy)
+ (ray.oz - sphere.cz) * (ray.oz - sphere.cz) - sphere.r2;
t1 = b * b - 4.0 * a * c;
if(t1 > 0)
{
t2 = sqrt(t1);
t3 = 2.0 * a;
close = -(b + t2) / t3;
farther = -(b - t2) / t3;
}
}
return (double)((close < farther) ? close : farther);
}
/*
** Return time to intersection with a plane
*/
double plane_intersect(struct RAY ray, struct PLANE plane)
{
double p1, p2, p3;
p1 = plane.px * plane.ny + plane.py * plane.ny + plane.pz * plane.nz;
p2 = ray.ox * plane.nx + ray.oy * plane.ny + ray.oz * plane.nz;
p3 = ray.dx * plane.nx + ray.dy * plane.ny + ray.dz * plane.nz;
return (double)((p1-p2)/p3);
}
/*
** Given a ray normal to the surface of reflection (n), an incident ray (i),
** return a reflected ray (r)
*/
void reflect(struct VECTOR *n, struct VECTOR *i, struct VECTOR *r)
{
double ndotn, idotn;
ndotn = (n->dx * n->dx + n->dy * n->dy + n->dz * n->dz);
idotn = (n->dx * i->dx + n->dy * i->dy + n->dz * i->dz);
r->dx = i->dx - (2.0 * (idotn) / ndotn) * n->dx;
r->dy = i->dy - (2.0 * (idotn) / ndotn) * n->dy;
r->dz = i->dz - (2.0 * (idotn) / ndotn) * n->dz;
}
int trace(void)
{
int x,y,c;
static struct PLANE plane = { 0.0, 1.0, 0.001, -6.0, 0.0, 0.0};
static struct SPHERE sphere = { 0.0, 2.0, 8.0, 49.0 };
struct RAY ray;
struct VECTOR v1, v2, v3;
double t1, t2, time;
c=0;
for(x=0;x<320 && !kbhit(); x+=2)
{
for(y = 0; y < 200; y+=2)
{
/* Ray needs to be reset because of
alterations in reflections */
ray.ox = 0.0;
ray.oy = 0.0;
ray.oz = 0.0;
ray.dz = 1.0;
ray.dy = -((double)y - 100.0) / 75.0;
ray.dx = ((double)x - 160.0) / 80.0;
t1 = sphere_intersect(ray,sphere);
t2 = plane_intersect(ray,plane);
if(t1 > 0.0 && (t2 < 0.0 || t2 > t1)) /* Circle in fore? */
{
v1.dx = ray.dx; v1.dy = ray.dy; v1.dz = ray.dz;
v2.dx = ((ray.dx * t1 + ray.ox) - sphere.cx);
v2.dy = ((ray.dy * t1 + ray.oy) - sphere.cy);
v2.dz = ((ray.dz * t1 + ray.oz) - sphere.cz);
reflect(&v2,&v1, &v3);
ray.ox += ray.dx * t1;
ray.oy += ray.dy * t1;
ray.oz += ray.dz * t1;
ray.dx = v3.dx;
ray.dy = v3.dy;
ray.dz = v3.dz;
t2 = plane_intersect(ray,plane);
if(t2 > 0.0)
{
int cr;
cr = (abs((int)(t2 * ray.dz + ray.oz)) % 2) +
2 * (abs((int)(t2 * ray.dx + ray.ox))
% 2);
plot(x,y,cr);
}
else plot(x,y,1);
}
else if(t2 > 0.0)
{
int cr;
cr = (abs((int)(t2 * ray.dz + ray.oz)) % 2) + 2 *
(abs((int)(t2 * ray.dx + ray.ox)) % 2);
plot(x,y,cr);
}
}
}
}
void main(void)
{
set_mode(4); /* Built for CGA, 4 color mode through BIOS */
trace();
getch();
set_mode(3);
}